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Abstract 

Using 99 4 field theory and Monte Carlo (MC) simulation we in- 
vestigate the finite-size effects of the magnetization M for the three- 
dimensional Ising model in a finite cubic geometry with periodic bound- 
ary conditions. The field theory with infinite cutoff gives a scal- 
ing form of the equation of state h/M s = f{hL^^ u ,t/h 1 ^ 5 ) where 
t = (T — T c )/T c is the reduced temperature, h is the external field 
and L is the size of system. Below T c and at T c the theory predicts 
a nonmonotonic dependence of f(x,y) with respect to x = hL^l v at 
fixed y = t/h 1 ^ 5 and a crossover from nonmonotonic to monotonic 
behaviour when y is further increased. These results are confirmed by 
MC simulation. The scaling function f(x,y) obtained from the field 
theory is in good quantitative agreement with the finite-size MC data. 
Good agreement is also found for the bulk value /(oo, 0) at T c . 
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1 Introduction 



The scaling equation of state near a critical point provides fundamental in- 
formation on the critical behaviour of a thermodynamic system. For bulk 
Ising-like systems accurate predictions have been made recently []J on the 
basis of the ip 4 field theory in three dimensions. Testing these predictions by 
Monte Carlo simulations would be of considerable interest. Such simu- 
lations, however, are necessarily made only for finite systems and thus phe- 
nomenological concepts like finite-size scaling |3j are needed to perform ex- 
trapolations from mesoscopic lattices. 

In order to test the theory in a more conclusive way it is desirable to go be- 
yond the phenomenological finite-size scaling concept and to calculate explic- 
itly the finite-size effects on the equation of state, i.e., on the magnetization 
as a function of the temperature T and external field h in a finite geometry. 
Although such field-theoretic calculations ||, |^, [F] are perturbative and not 
exact they have been found to be in good agreement with the MC simulations. 
So far the calculations of thermodynamic quantities were restricted to the 
case of zero external field h. At finite h, only the order-parameter distribu- 
tion function was calculated || [| |Hj. In the present paper we extend these 
calculations to the equation of state of the three-dimensional Ising model 
in a finite cubic geometry at finite h. For simplicity these calculations are 
performed at infinite cutoff and thus our results neglect lattice effects. The 
latter have been shown [[0], |l2j to yield (exponentially small) non-universal 



non-scaling contributions. Here we focus our interest on the universal scal- 
ing part of the equation of state. Our theory predicts non-monotonic effects 
in the h dependence of the equation of state which we then confirm with 
surprisingly good agreement by standard Monte Carlo simulations. 

The identification of such non-monotonic effects is of great practical impor- 
tance. If, for example, the critical temperature T C (L) of a lattice with L 3 
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sites in three dimensions varies asymptotically as T C (L) — T c (oo) oc 1/L y 
with some correlation length exponent y — then a plot of the numerical 
T C (L) versus 1/L y gives the extrapolated T c (oo) as an intercept. If, however, 
higher order terms make the curve T C (L) non-monotonic, then such a plot for 
finite L, where the non-monotonicity is not yet visible, would give a wrong 
estimate for T c (oo). Similar effects would make estimates of other quantities 
unreliable, and there exist examples of such estimates in the literature. 



2 Field-theoretic calculations 

We consider the <£> 4 model with the standard Landau- Ginzburg- Wilson Hamil- 
tonian 



H(h) 



v 



\ r W 2 + t;(w) 2 + ?W 4 -h-<p 



(1) 



where y(x) is a one-component field in a finite cube of volume V = L d and 
h is a homogeneous external field. We assume periodic boundary conditions. 
Accordingly we have 

y?(x) =L- d $> k e ik - x (2) 



k 

2tt, 



where the summation J2k runs over discrete k = vectors with componets 
kj = jr-m,j,mj = 0, ±1, ±2, j = 1, 2, d in the range —A < kj < A, i.e., 
with a sharp cutoff A. The temperature enters through r = r 0c + dot, i = 
(T-T c )/T c . 



As pointed out recently [|T^] the Hamiltonian (1) for periodic boundary con- 
ditions with a sharp cutoff A does not correctly describe the exponential 
size dependence of physical quantities of finite lattice models in the region 
£ ^> L. Instead of (1), a modified continuum Hamiltonian with a smooth 
cutoff (T^j would be more appropriate. Even better would be to employ the 
lattice version of the f A theory to describe the non-scaling lattice effects of 
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finite Ising models in the region £ ^> L. In the present paper, however, we 
shall neglect such effects by taking the limit A — > oo (see below). 

The fluctuating homogeneous part of the order-parameter of the Hamiltonian 
(1) is $ = V~ l j v d d xip(x) = L~ d (p Q . As previously [5] (p is decomposed as 

y?(x) = $ + o-(x) (3) 

where <r(x) includes all inhomogeneous modes 

*M = L- d ^ ^ k e lk - x (4) 

The order-parameter distribution function P($) = P(3>, t, /i, L) is defined by 
functional integration over a, 

P($, t, h, L) = Z(h)- 1 J Dae~ H{h) , (5) 



where 



Z(h) = f°° I Dae- H{h) (6) 

J —oo J 



is the partition function of system. This distribution function depends also 
on the cutoff A which implies non-scaling finite-size effects [JTTJ, [12| . From the 
order-parameter distribution function P(<3?) we can calculate the magnetiza- 
tion 

("OO 

M=<|$|>=/ d$|$|P($). (7) 



— oo 



The functional integration over a in Eqs. (5) and (6) can only be done per- 
turbatively (|, ||, §]. Recently a novel perturbation approach was presented 
[0, ||. Using this approach the order-parameter distribution can be written 
in the form 

P($) = e -^ e// W / /°° d$e" He// ( # ) (8) 

J — oo 
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where the (bare) effective Hamiltonian of the Ising-like system reads 

»'"(*) = W)-;EM#|f}. 0» 

#o(<M) = L rf (^ro$ 2 + Wo$ 4 -/i-$), (10) 

with r 0L = r + 12u $ 2 and y m(r) = (2L d /3M ) 1/2 (r + frm 2 ). The function 
in Eq.(9) is defined as 

r°° 1 
Z 1 [y}= ds sexp(--ys 2 - s 4 ). (11) 
jo 2 



After minimal renormalization at fixed dimension d < 4 JIB) in the limit A 



oo, we obtain the finite-size scaling form of the order-parameter distribution 
function g, g] 

P($, t, h, L) = L l3/ "p{hL /35/ ", tL 1/u , $1//"), (12) 

where 

/ \ exp[-F(z,g,z)] 
p{x, q, z) = — — r — -, (13) 

with x = hL ps ' v , q = tL 1 ^, z = and 

F(x,q,z) = C2(x, q)z 2 + c±(x, q)z 4 — x ■ z 
_1 l r Z 1 [y m (f L (x,q,z))] 

Here g = Q*t(L/^ ) 1/i/ and z = (2Q*) /3 ($/A M )(^/^ ) /3/l/ are dimensionless 
variables normalized by the asymptotic amplitudes £o and Am of the bulk cor- 
relation length £ = £ s ot~ u at h — above T c and of the bulk order-parameter 
Mbuik = Am^I 13 at /i = below T c . The bulk parameter Q* is known . 
The coefficients c 2 (x, q) and c 4 (x, g) read for d = 3 

c 2 (x,g) = (647TM*)- 1 g£(x,g) 3 -( 2/3+1 ^(l + 12m*), (15) 
c 4 (x,g) = {256nu*)- 1 £(x,q) 3 ~ 4/3/u {l + 36u*), (16) 
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where u* is the fixed point value of the renormalized coupling |Tj|. In three 
dimensions we have 



y m (r L (x,q,z)) = [67ru*£(x,q)}- l ^ L (x,q,z)£(x,q) 2 + An 2 in 2 ], (17) 
r L (x,q,z) = ql(xA)~ 1/v + (3/2)£(x,q)~ 2 ^z 2 . (18) 

The auxiliary scaling function £(x, q) of the flow parameter is determined by 

£(xA?' 2 = (^u*) 1/2 [y(x,q) + 12Uy(x,q),x)], (19) 
y(x,q) = (Anu*)- 1/2 ((x,q) 3/2 - 1/v q, (20) 
x = A M (2QT^o /u (^u*) 1 ^V8£(x,q)^- 3 ^x, (21) 

where 

$jy £ \ = ^0°° ds g2 cosh (^) exp(-|yg 2 - s A ) ^ 
J °° ds cosh(xs) exp(— \ys 2 — s 4 ) 

From Eqs. (7) and (12) we obtain the scaling form 

M(h,t,L) = L-^f M (hL^\tL^), (23) 

f M (x,q) = S-oodz\z\exp[-F(^Q^)] ^ (24) 
J^o(izexp[-F(x,g,z)] 

Correspondingly the asymptotic equation of state for a finite Ising-like system 
in the limit of zero lattice spacing [[12] can be written as 

h/M s = f(hL p5/ \ t/h l/ps ) (25) 

where 

f(x,y) = x/f M (x,yx 1 ^ 5 ) 5 . (26) 

In the comparison with the MC data of the simple-cubic (sc) Ising model in 
Section 4, the quantities h, M and L are used in a dimensionless form in units 
of the lattice constant, see Sect. 4 of Ref. 0. We have taken the bulk param- 



eters u* = 0.0412, Q* = 0.945 from Ref.H. From Ref.pi] we have taken the 
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bulk amplitudes £o = 0.495, Am = 1-71 in units of the lattice constant (of 
the sc Ising model) and the bulk critical exponents (3 = 0.3305, v = 0.6335. 
Thus our determination of the scaling function f(x, y) does not require a new 
adjustment of nonuniversal parameters. 

Taking the limit hL /35 / u — > oo at fixed t/h 1 ^ 5 , we obtain the scaling form of 
the bulk equation of state 

h/M 5 = hit/h 1 ^ 8 ) = /(oo, t/h 1 ^ 5 ). (27) 

At T = T c we find from Eqs. (23)- (27) 

h/M 5 = D C = /(oo, 0) = 0.202 (28) 

in three dimensions (in units of the lattice constant). 



3 Monte Carlo simulation 

Standard heat bath techniques were used for the Glauber kinetic Ising model, 
with multi-spin coding (16 spins in each 64-bit computer word). Since the 
effect could be seen best in lattices of intermediate size L ~ 80 for L x L x 
L spins, memory requirements were tiny and only trivial parallelization by 
replication, not by domain decomposition, was used. However, thousands of 
hours of processor time were needed since we see the non-monotonic effects 
clearly in our figure if h/M 5 is plotted. The exponent 5 is nearly five and thus 
five percent accuracy in h/M 5 requires one percent accuracy in the directly 
simulated magnetization M. 



Earlier simulations by the same author and algorithm [fL5[] did not show the 
non-monotonic effects at T = T c because they were not searched for; at that 
time the field-theoretical predictions presented above were not yet known. 
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But even if they had been known it is doubtful that with the Intel Paragon 
used in |L5| instead of the Cray-T3E now the non-monotonic trends would 
have been seen in about the same computer time. 



Because of the limited system size, errors in the magnetization of order 1 
+ const /L 13 ^ ~ 1 + const/\/L are expected; in light of these errors the 
agreement to be presented now is surprisingly good. 



4 Results and discussion 

Figure 1 compares our Monte Carlo results for 32 3 and 80 3 spins with our 
theoretical predictions, Eqs. (23)-(26), and shows good quantitative agree- 
ment; below T c a pronounced peak is seen in the scaling function (Fig. la), at 
T c it is somewhat weaker (Fig. lb), above T c it is barely visible (Fig.lc), and 
far above T c (Fig. Id) it has vanished. The simulations, which at T = T c have 
an accuracy of the order of one and five percent for L = 32 and 80, respec- 
tively, agree nicely with the theoretical predictions. (Far above T c therefore 
no simulations were made.) 

Thus, if one tries to determine the bulk critical amplitude of h/M s at the 
critical isotherm, then with varying L at a fixed field one first gets a too small 
value (left border of the figures), then a too high value (peak), and then a 
roughly correct value (plateau in the right part of the figure). 

These non-monotonicities do not vanish if we take the lattices large enough. 
They are part of the asymptotic scaling function and thus whatever lattice 
size L we take there will be a field h oc l/L 135 ^ where h/M 5 has a maximum 
and near which it thus varies non-monotonically with L. 



We also tested the prediction of [11, 12] that for fixed T < T c the leading 
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finite-size deviation from the bulk value of M should vanish exponentially 
in L, and not with a power law oc 1/L (as predicted by perturbation the- 
ory based on the separation of the zero- mode [4-7]). Fig. 2a shows again 
non-monotonic behaviour, in both three and five dimensions. But only in 
three dimensions these data are accurate enough to distinguish between a 
tail varying exponentially and one cx 1/L d ; Fig.2b clearly supports the the- 



oretically predicted [11, 12 1 exponential variation. (These three-dimensional 



data were taken with Ito's fast algorithm 



In order to make contact with bulk properties we also used simulations at 
the critical isotherm with 1292 3 spins ||15||. From the simulations we obtain 



the bulk value of h/M 5 as D c = 0.21 ± 0.02. Our field-theoretic result in Eq. 
(28) is in very good agreement with this value. 

It is interesting to compare our simulation result also with other bulk theories. 



From series expansion Zinn and Fisher |L7| obtained C c = 0.299 for the 
amplitude of the bulk susceptibility at the critical isotherm \ = C c \h\~' y ^ s 
with 7 = 1.2395, v = 0.6320. This leads to D c = 0.182 according to the 
relation D c = {C c 5)~ & . 



D c is also contained in the universal combination of amplitudes [[16 

R x = YD c A 5 m 1 (29) 

where Y is the amplitude of the bulk susceptibility x — T£~ 7 above T c at h = 
and Am is the amplitude of the spontaneous magnetization M^uik = AmI^ 
below T c . Using (p 4 field theory at d — 3 dimensions Guida and Zinn- Justin 
[H] have obtained R x = 1.649. They also used the e = 4 — d expansion 
and obtained R x = 1.674 at e = 1. Using these values for R x and the 
high-temperature series expansion results |14[ T = 1.0928, Am = 1-71 and 



S = (dp + i)/{dv - 7) with 7 = 1.2395, v = 0.6335, we obtain from Eq. (29) 
D c = 0.205 (d = 3 field theory) and D c = 0.202 (e-expansion), respectively, 
in good agreement with our theoretical result, Eq. (28), and with our MC 
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simulation. 



In summary, our simulations confirmed a posteriori our theoretical predic- 
tions of Sect. 2 for the asymptotic finite-size effects in the three-dimensional 
Ising model. The agreement in Fig.f is remarkable in view of the fact that the 
non-universal parameters of the theory were adjusted only to bulk parameters 
of the Ising model at h — and not to any finite-size MC data. 
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Figure Captions 



Fig. 1. Scaling plot (in units of the lattice constant, see Ref. ||) of h/M s 
versus x = hL ps/v below T c with t/h 1/f35 = -1.0 in part a, t/h 1 ^ 5 = (that 
means T = T c ) in part b, t/h 1 ^ 5 = 1.0 in part c, and t/h 1 ^ 6 = 1.6 in part d. 
Monte Carlo data for L = 32 and 80. Solid line is the theoretical prediction 
Eqs. (23)- (26); no Monte Carlo data are shown in part d, where h/M s is a 
monotonic function of external field x = hL^l v . 

Fig. 2. a) Monte Carlo data for the spontaneous magnetization in units 
of the lattice constant in three (diamonds and plusses) and five (square) 
dimensions at T/T c = 0.99, versus linear lattice size L. The horizontal line 
M = 0.3671 for three dimensions is determined from L = 2496. 



b) Selected three-dimensional data from part a) are shown as M(L = oo) — 
M(L) versus L in a semilogarithmic plot. The straight line represents an 
exponential decay |TT[ [12], the two curved lines are power law decays 1/L 2 
and 1/L 3 which fail to fit the data. The five-dimensional data of part a) 
were not accurate enough to distinguish between an exponential and a 1/L 5 
decay. 
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